setwd("~/Documents/GitHub/Resultados/docs/PrimeroDiffExpAllResults/Clasificando/Abundances")
load('CheackPointOne.RData')

Quinta Clasificacion: pEhExvsUmasM

head(pEhExvsUmasM,10);
##          GenId   UmasM_1    UmasM_2    UmasM_3    pEhEx_1    pEhEx_2    pEhEx_3
## 1  EHI_000130A  52.79569  428.31302  877.32975  175.05490  432.49403  478.29632
## 2  EHI_000140A 117.59041  370.70763   43.71038  353.07682  340.82410   38.78078
## 3  EHI_000240A 743.93931 1164.98430 7104.49766 1174.94473 1932.12006 2122.17063
## 4  EHI_000250A 704.34254  212.80109 1578.25690  324.89002  326.72103  658.55515
## 5  EHI_000260A 118.79031   96.23489  138.93656  109.03843   66.98956   50.27139
## 6  EHI_000280A  57.59530   59.63852   34.34387   54.89009   83.44314   61.76199
## 7  EHI_000290A  26.39785   19.65360   70.24882   18.54395   14.10307   81.15238
## 8  EHI_000300A  80.39344  115.21078   15.61085  133.51645  115.17504   31.59916
## 9  EHI_000410A  20.39834   24.39758  138.93656   28.18681   24.68037   65.35280
## 10 EHI_000430A  31.19746   35.91866   14.04976   25.96153   14.10307   10.77244
nbreaks <- 10
data5 <- pEhExvsUmasM;       head(data5)
##         GenId   UmasM_1    UmasM_2    UmasM_3    pEhEx_1    pEhEx_2    pEhEx_3
## 1 EHI_000130A  52.79569  428.31302  877.32975  175.05490  432.49403  478.29632
## 2 EHI_000140A 117.59041  370.70763   43.71038  353.07682  340.82410   38.78078
## 3 EHI_000240A 743.93931 1164.98430 7104.49766 1174.94473 1932.12006 2122.17063
## 4 EHI_000250A 704.34254  212.80109 1578.25690  324.89002  326.72103  658.55515
## 5 EHI_000260A 118.79031   96.23489  138.93656  109.03843   66.98956   50.27139
## 6 EHI_000280A  57.59530   59.63852   34.34387   54.89009   83.44314   61.76199

Log-NormalizaciĂ³n

sample1   <- data5$pEhEx_1; sample2   <- data5$pEhEx_2; sample3   <- data5$pEhEx_3;
samplevs1 <- data5$UmasM_1;  samplevs2 <- data5$UmasM_2;  samplevs3 <- data5$UmasM_3;
log2sample1 <- log2(sample1+1);         log2sample2 <- log2(sample2+1)
log2sample3 <- log2(sample3+1);         log2samplevsumasM1 <- log2(samplevs1+1)
log2samplevsumasM2 <- log2(samplevs2+1); log2samplevsumasM3 <- log2(samplevs3+1)
data5 <- cbind(data5, log2sample1,log2sample2,log2sample3,
               log2samplevsumasM1,log2samplevsumasM2,log2samplevsumasM3)
head(data5)
##         GenId   UmasM_1    UmasM_2    UmasM_3    pEhEx_1    pEhEx_2    pEhEx_3
## 1 EHI_000130A  52.79569  428.31302  877.32975  175.05490  432.49403  478.29632
## 2 EHI_000140A 117.59041  370.70763   43.71038  353.07682  340.82410   38.78078
## 3 EHI_000240A 743.93931 1164.98430 7104.49766 1174.94473 1932.12006 2122.17063
## 4 EHI_000250A 704.34254  212.80109 1578.25690  324.89002  326.72103  658.55515
## 5 EHI_000260A 118.79031   96.23489  138.93656  109.03843   66.98956   50.27139
## 6 EHI_000280A  57.59530   59.63852   34.34387   54.89009   83.44314   61.76199
##   log2sample1 log2sample2 log2sample3 log2samplevsumasM1 log2samplevsumasM2
## 1    7.459882    8.759868    8.904774           5.749419           8.745886
## 2    8.467919    8.417110    5.314000           6.889844           8.538024
## 3   10.199605   10.916716   11.052005           9.540979          10.187333
## 4    8.348241    8.356324    9.365349           9.462180           7.740125
## 5    6.781864    6.087241    5.680082           6.904367           6.603402
## 6    5.804521    6.399908    5.971819           5.872713           5.922163
##   log2samplevsumasM3
## 1           9.778619
## 2           5.482538
## 3          12.794720
## 4          10.625030
## 5           7.128629
## 6           5.143388
save.image('CheckPointFifth.RData')
setwd("~/Documents/GitHub/Resultados/docs/PrimeroDiffExpAllResults/Clasificando/Abundances")
#load('CheckPointTwo.RData')
library(ggplot2);library(dplyr);library("fitdistrplus");
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
## Loading required package: MASS
## 
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
## 
##     select
## Loading required package: survival
library("MASS");library("survival")
head(data5)
##         GenId   UmasM_1    UmasM_2    UmasM_3    pEhEx_1    pEhEx_2    pEhEx_3
## 1 EHI_000130A  52.79569  428.31302  877.32975  175.05490  432.49403  478.29632
## 2 EHI_000140A 117.59041  370.70763   43.71038  353.07682  340.82410   38.78078
## 3 EHI_000240A 743.93931 1164.98430 7104.49766 1174.94473 1932.12006 2122.17063
## 4 EHI_000250A 704.34254  212.80109 1578.25690  324.89002  326.72103  658.55515
## 5 EHI_000260A 118.79031   96.23489  138.93656  109.03843   66.98956   50.27139
## 6 EHI_000280A  57.59530   59.63852   34.34387   54.89009   83.44314   61.76199
##   log2sample1 log2sample2 log2sample3 log2samplevsumasM1 log2samplevsumasM2
## 1    7.459882    8.759868    8.904774           5.749419           8.745886
## 2    8.467919    8.417110    5.314000           6.889844           8.538024
## 3   10.199605   10.916716   11.052005           9.540979          10.187333
## 4    8.348241    8.356324    9.365349           9.462180           7.740125
## 5    6.781864    6.087241    5.680082           6.904367           6.603402
## 6    5.804521    6.399908    5.971819           5.872713           5.922163
##   log2samplevsumasM3
## 1           9.778619
## 2           5.482538
## 3          12.794720
## 4          10.625030
## 5           7.128629
## 6           5.143388

Muestra 1

log2sample1 <- data5$log2sample1; head(mean(log2sample1)); head(sd(log2sample1))
## [1] 6.24652
## [1] 2.884429
head(log2sample1,5)
## [1]  7.459882  8.467919 10.199605  8.348241  6.781864
summary(data5)
##     GenId              UmasM_1            UmasM_2             UmasM_3         
##  Length:4919        Min.   :     0.0   Min.   :     0.00   Min.   :      0.0  
##  Class :character   1st Qu.:    16.8   1st Qu.:    18.30   1st Qu.:     14.0  
##  Mode  :character   Median :    46.8   Median :    52.18   Median :     54.6  
##                     Mean   :  1915.1   Mean   :  1003.67   Mean   :   3444.8  
##                     3rd Qu.:   198.0   3rd Qu.:   193.15   3rd Qu.:    295.0  
##                     Max.   :340994.2   Max.   :145896.83   Max.   :1488475.8  
##     pEhEx_1             pEhEx_2             pEhEx_3          log2sample1    
##  Min.   :     0.00   Min.   :     0.00   Min.   :     0.0   Min.   : 0.000  
##  1st Qu.:    16.32   1st Qu.:    15.28   1st Qu.:    15.1   1st Qu.: 4.114  
##  Median :    48.96   Median :    48.19   Median :    52.4   Median : 5.643  
##  Mean   :  1394.44   Mean   :  1752.18   Mean   :  1898.2   Mean   : 6.247  
##  3rd Qu.:   199.53   3rd Qu.:   222.71   3rd Qu.:   231.2   3rd Qu.: 7.648  
##  Max.   :213518.02   Max.   :279409.95   Max.   :724577.3   Max.   :17.704  
##   log2sample2      log2sample3     log2samplevsumasM1 log2samplevsumasM2
##  Min.   : 0.000   Min.   : 0.000   Min.   : 0.000     Min.   : 0.000    
##  1st Qu.: 4.025   1st Qu.: 4.007   1st Qu.: 4.154     1st Qu.: 4.270    
##  Median : 5.620   Median : 5.739   Median : 5.579     Median : 5.733    
##  Mean   : 6.187   Mean   : 6.152   Mean   : 6.222     Mean   : 6.256    
##  3rd Qu.: 7.805   3rd Qu.: 7.860   3rd Qu.: 7.637     3rd Qu.: 7.601    
##  Max.   :18.092   Max.   :19.467   Max.   :18.379     Max.   :17.155    
##  log2samplevsumasM3
##  Min.   : 0.000    
##  1st Qu.: 3.912    
##  Median : 5.798    
##  Mean   : 6.317    
##  3rd Qu.: 8.210    
##  Max.   :20.505
ndata5    <- length(log2sample1)
hist(log2sample1, breaks = nbreaks, col= rainbow(25,0.3), 
     main = 'Log2 sample1')

meanlog2sample1 <- mean(log2sample1); head(meanlog2sample1)
## [1] 6.24652
StdDevlog2sample1 <- sd(log2sample1); head(StdDevlog2sample1)
## [1] 2.884429
Normlog2sample1 <- (log2sample1-meanlog2sample1)/StdDevlog2sample1; head(Normlog2sample1)
## [1]  0.4206592  0.7701346  1.3704911  0.7286438  0.1855978 -0.1532363
tst<- Normlog2sample1

hist(tst, breaks = nbreaks, col= 1:5, 
     main = 'Normalized Log2sample1',
     xlab='pEhEx1',
     ylab= 'Frequency pEhEx')

fw1<-fitdist(tst, "norm")
plotdist(tst, histo = TRUE, demp = TRUE)

nnorm.f <- fitdist(tst,"norm")
summary(nnorm.f)
## Fitting of the distribution ' norm ' by maximum likelihood 
## Parameters : 
##           estimate Std. Error
## mean -1.310539e-16 0.01425665
## sd    9.998983e-01 0.01008093
## Loglikelihood:  -6979.259   AIC:  13962.52   BIC:  13975.52 
## Correlation matrix:
##      mean sd
## mean    1  0
## sd      0  1
par(mfrow=c(2,2))
denscomp(nnorm.f,legendtext = 'Dist Normal')
qqcomp(nnorm.f,legendtext = 'Dist Normal')
cdfcomp(nnorm.f,legendtext = 'Dist Normal')
ppcomp(nnorm.f,legendtext = 'Dist Normal')

probs <- c();
probs[8] = 0.175;  probs[9] = 0.825; 
probs[7] = 0.15;   probs[10] = 0.85;   
probs[6] = 0.125;  probs[11] = 0.875; 
probs[5] = 0.1;    probs[12] = 0.9;    
probs[4] = 0.075;  probs[13] = 0.925; 
probs[3] = 0.05;   probs[14] = 0.95;   
probs[2] = 0.025;  probs[15] = 0.975; 
probs[1] = 0.005;  probs[16] = 0.995;  
CuantilesData <- quantile(tst,prob = probs)
CuantilesModel <- qnorm(probs, mean=0, sd=1)
Cuantilillos <- t(CuantilesModel)
colnames(Cuantilillos) <- c('0.5%','2.5%','5%','7.5%',
                            '10%','12.5%','15%','17.5%',
                            '82.5%','85%','87.5%','90%',
                            '92.5%','95%','97.5%','99.5%')
Cuantilillos <- t(Cuantilillos)
colnames(Cuantilillos) <- c('Cuantiles Ajuste')
print(Cuantilillos)
##       Cuantiles Ajuste
## 0.5%        -2.5758293
## 2.5%        -1.9599640
## 5%          -1.6448536
## 7.5%        -1.4395315
## 10%         -1.2815516
## 12.5%       -1.1503494
## 15%         -1.0364334
## 17.5%       -0.9345893
## 82.5%        0.9345893
## 85%          1.0364334
## 87.5%        1.1503494
## 90%          1.2815516
## 92.5%        1.4395315
## 95%          1.6448536
## 97.5%        1.9599640
## 99.5%        2.5758293

CĂ¡lculo de cuantiles

CuantilesA <- matrix(0,8,2)
CuantilesD <- matrix(0,8,2)
colnames(CuantilesA) <- c('LimInf','LimSup')
colnames(CuantilesD) <- c('LimInf','LimSup')
rownames(CuantilesA) <- c('65','70','75','80','85','90','95','99')
rownames(CuantilesD) <- c('65','70','75','80','85','90','95','99')
CuantilesA[1,1] <-CuantilesData[8]; CuantilesA[1,2] <-CuantilesData[9]
CuantilesA[2,1] <-CuantilesData[7]; CuantilesA[2,2] <-CuantilesData[10]
CuantilesA[3,1] <-CuantilesData[6]; CuantilesA[3,2] <-CuantilesData[11]
CuantilesA[4,1] <-CuantilesData[5]; CuantilesA[4,2] <-CuantilesData[12]
CuantilesA[5,1] <-CuantilesData[4]; CuantilesA[5,2] <-CuantilesData[13]
CuantilesA[6,1] <-CuantilesData[3]; CuantilesA[6,2] <-CuantilesData[14]
CuantilesA[7,1] <-CuantilesData[2]; CuantilesA[7,2] <-CuantilesData[15]
CuantilesA[8,1] <-CuantilesData[1]; CuantilesA[8,2] <-CuantilesData[16]
CuantilesD[1,1] <-Cuantilillos[8];  CuantilesD[1,2] <-Cuantilillos[9]
CuantilesD[2,1] <-Cuantilillos[7];  CuantilesD[2,2] <-Cuantilillos[10]
CuantilesD[3,1] <-Cuantilillos[6];  CuantilesD[3,2] <-Cuantilillos[11]
CuantilesD[4,1] <-Cuantilillos[5];  CuantilesD[4,2] <-Cuantilillos[12]
CuantilesD[5,1] <-Cuantilillos[4];  CuantilesD[5,2] <-Cuantilillos[13]
CuantilesD[6,1] <-Cuantilillos[3];  CuantilesD[6,2] <-Cuantilillos[14]
CuantilesD[7,1] <-Cuantilillos[2];  CuantilesD[7,2] <-Cuantilillos[15]
CuantilesD[8,1] <-Cuantilillos[1];  CuantilesD[8,2] <-Cuantilillos[16]
print(CuantilesD)
##        LimInf    LimSup
## 65 -0.9345893 0.9345893
## 70 -1.0364334 1.0364334
## 75 -1.1503494 1.1503494
## 80 -1.2815516 1.2815516
## 85 -1.4395315 1.4395315
## 90 -1.6448536 1.6448536
## 95 -1.9599640 1.9599640
## 99 -2.5758293 2.5758293
print(CuantilesA)
##        LimInf    LimSup
## 65 -0.8877976 0.8234548
## 70 -0.9174931 0.9668385
## 75 -0.9827617 1.1346553
## 80 -1.0188953 1.3615138
## 85 -1.1000839 1.6923518
## 90 -1.1462228 2.1593960
## 95 -1.3906263 2.6069712
## 99 -2.1655999 3.1547856
col_sequence <- rainbow(n = 7, alpha = 0.35, start = 0, end = 1)
par(mfrow=c(2,1))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2 pEhEx - DATA', lty = 9)

abline(v=CuantilesA[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesA[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesA[2,1], lty=2, col="darkblue"); # 70% INFERIOR
abline(v=CuantilesA[2,2], lty=2, col="darkblue") # 70% SUPERIOR
abline(v=CuantilesA[3,1], lty=2, col="aquamarine4"); # 75% INFERIOR
abline(v=CuantilesA[3,2], lty=2, col="aquamarine4"); # 75% SUPERIOR
abline(v=CuantilesA[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesA[4,2], lty=2, col="green"); # 80% SUPERIOR
abline(v=CuantilesA[5,1], lty=2, col="brown"); # 85% INFERIOR
abline(v=CuantilesA[5,2], lty=2, col="brown"); # 85% SUPERIOR
abline(v=CuantilesA[6,1], lty=2, col="red");  # 90% INFERIOR
abline(v=CuantilesA[6,2], lty=2, col="red"); # 90% SUPERIOR
abline(v=CuantilesA[7,1], lty=2, col="blue");  # 95% INFERIOR
abline(v=CuantilesA[7,2], lty=2, col="blue"); # 95% SUPERIOR
abline(v=CuantilesA[8,1], lty=2, col="orange");  # 99% INFERIOR
abline(v=CuantilesA[8,2], lty=2, col="orange"); # 99% SUPERIOR
legend("topright",
       legend=c("65%","70%","75%","80%","85%","90%","95%","99%"), 
       pch=c(1,2,3,4,5,6,7,8),
       col=c("darkgoldenrod4","darkblue","aquamarine4",
             "green", "brown","red","blue","orange"))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2   pEhEx - ADJUSTED', lty=9)
abline(v=CuantilesD[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesD[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesD[2,1], lty=2, col="darkblue");  # 70% INFERIOR
abline(v=CuantilesD[2,2], lty=2, col="darkblue"); # 70% SUPERIOR
abline(v=CuantilesD[3,1], lty=2, col="aquamarine4");  # 75% INFERIOR
abline(v=CuantilesD[3,2], lty=2, col="aquamarine4"); # 75% SUPERIOR
abline(v=CuantilesD[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesD[4,2], lty=2, col="green"); # 80% SUPERIOR
abline(v=CuantilesD[5,1], lty=2, col="brown");  # 85% INFERIOR
abline(v=CuantilesD[5,2], lty=2, col="brown"); # 85% SUPERIOR
abline(v=CuantilesD[6,1], lty=2, col="red");  # 90% INFERIOR
abline(v=CuantilesD[6,2], lty=2, col="red"); # 90% SUPERIOR
abline(v=CuantilesD[7,1], lty=2, col="blue");  # 95% INFERIOR
abline(v=CuantilesD[7,2], lty=2, col="blue"); # 95% SUPERIOR
abline(v=CuantilesD[8,1], lty=2, col="orange");  # 99% INFERIOR
abline(v=CuantilesD[8,2], lty=2, col="orange"); # 99% SUPERIOR
legend("topright",
       legend=c("65%","70%","75%","80%","85%","90%","95%","99%"), 
       pch=c(1,2,3,4,5,6,7,8),
       col=c("darkgoldenrod4","darkblue","aquamarine4",
             "green", "brown","red","blue","orange"))

Grafica Cuantiles del \(65\%\) y \(80\%\)

par(mfrow=c(2,1))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2  pEhEx - DATA', lty=9)

abline(v=CuantilesA[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesA[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesA[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesA[4,2], lty=2, col="green"); # 80% SUPERIOR
legend("topright",legend=c("65%","80%"), 
       pch=c(1,2),col=c("darkgoldenrod4","green"))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2  BaseMean - ADJUSTED', lty=9)

abline(v=CuantilesD[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesD[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesD[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesD[4,2], lty=2, col="green"); # 80% SUPERIOR
legend("topright",legend=c("65%","80%"), 
       pch=c(1,2),col=c("darkgoldenrod4","green"))

Grafica Cuantiles del \(70\%\) y \(85\%\)

par(mfrow=c(2,1))
hist(tst, breaks = nbreaks, col = col_sequence, 
     main = 'Normalized Log2 pEhEx - DATA', lty=9)

abline(v=CuantilesA[2,1], lty=2, col="darkblue"); 
abline(v=CuantilesA[2,2], lty=2, col="darkblue")
abline(v=CuantilesA[5,1], lty=2, col="brown"); 
abline(v=CuantilesA[5,2], lty=2, col="brown")
legend("topright",legend=c("70%","85%"),
       pch=c(1,2),#3,4,5,6,7,8),
       col=c("brown"))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2  pEhEx - ADJUSTED', lty=9)

abline(v=CuantilesD[2,1], lty=2, col="darkblue"); 
abline(v=CuantilesD[2,2], lty=2, col="darkblue")
abline(v=CuantilesD[5,1], lty=2, col="brown"); 
abline(v=CuantilesD[5,2], lty=2, col="brown")
legend("topright",legend=c("70%","85%"),
       pch=c(1,2),#3,4,5,6,7,8),
       col=c("brown"))

Muestra 2

log2sample2 <- data5$log2sample2; head(mean(log2sample2)); head(sd(log2sample2))
## [1] 6.18657
## [1] 3.14197
head(log2sample2,5)
## [1]  8.759868  8.417110 10.916716  8.356324  6.087241
ndata5    <- length(log2sample2)
hist(log2sample2, breaks = nbreaks, col= rainbow(25,0.3), 
     main = 'Log2 sample2')

Log-normalizacion

meanlog2sample2 <- mean(log2sample2); head(meanlog2sample2)
## [1] 6.18657
StdDevlog2sample2 <- sd(log2sample2); head(StdDevlog2sample2)
## [1] 3.14197
Normlog2sample2 <- (log2sample2-meanlog2sample2)/StdDevlog2sample2; head(Normlog2sample2)
## [1]  0.81900802  0.70991786  1.50547137  0.69057145 -0.03161338  0.06789964
tst<- Normlog2sample2
hist(tst, breaks = nbreaks, col= 1:5, 
     main = 'Normalized Log2sample1',
     xlab='pEhEx1',
     ylab= 'Frequency pEhEx')

Ajuste de modelo

fw1<-fitdist(tst, "norm")
plotdist(tst, histo = TRUE, demp = TRUE)

nnorm.f <- fitdist(tst,"norm")
summary(nnorm.f)
## Fitting of the distribution ' norm ' by maximum likelihood 
## Parameters : 
##           estimate Std. Error
## mean -1.093861e-16 0.01425665
## sd    9.998983e-01 0.01008093
## Loglikelihood:  -6979.259   AIC:  13962.52   BIC:  13975.52 
## Correlation matrix:
##      mean sd
## mean    1  0
## sd      0  1
par(mfrow=c(2,2))
denscomp(nnorm.f,legendtext = 'Dist Normal')
qqcomp(nnorm.f,legendtext = 'Dist Normal')
cdfcomp(nnorm.f,legendtext = 'Dist Normal')
ppcomp(nnorm.f,legendtext = 'Dist Normal')

probs <- c();
probs[8] = 0.175;  probs[9] = 0.825; 
probs[7] = 0.15;   probs[10] = 0.85;   
probs[6] = 0.125;  probs[11] = 0.875; 
probs[5] = 0.1;    probs[12] = 0.9;    
probs[4] = 0.075;  probs[13] = 0.925; 
probs[3] = 0.05;   probs[14] = 0.95;   
probs[2] = 0.025;  probs[15] = 0.975; 
probs[1] = 0.005;  probs[16] = 0.995;  
CuantilesData <- quantile(tst,prob = probs)
CuantilesModel <- qnorm(probs, mean=0, sd=1)
Cuantilillos <- t(CuantilesModel)
colnames(Cuantilillos) <- c('0.5%','2.5%','5%','7.5%',
                            '10%','12.5%','15%','17.5%',
                            '82.5%','85%','87.5%','90%',
                            '92.5%','95%','97.5%','99.5%')
Cuantilillos <- t(Cuantilillos)
colnames(Cuantilillos) <- c('Cuantiles Ajuste')
print(Cuantilillos)
##       Cuantiles Ajuste
## 0.5%        -2.5758293
## 2.5%        -1.9599640
## 5%          -1.6448536
## 7.5%        -1.4395315
## 10%         -1.2815516
## 12.5%       -1.1503494
## 15%         -1.0364334
## 17.5%       -0.9345893
## 82.5%        0.9345893
## 85%          1.0364334
## 87.5%        1.1503494
## 90%          1.2815516
## 92.5%        1.4395315
## 95%          1.6448536
## 97.5%        1.9599640
## 99.5%        2.5758293
CuantilesA <- matrix(0,8,2)
CuantilesD <- matrix(0,8,2)
colnames(CuantilesA) <- c('LimInf','LimSup')
colnames(CuantilesD) <- c('LimInf','LimSup')
rownames(CuantilesA) <- c('65','70','75','80','85','90','95','99')
rownames(CuantilesD) <- c('65','70','75','80','85','90','95','99')
CuantilesA[1,1] <-CuantilesData[8]; CuantilesA[1,2] <-CuantilesData[9]
CuantilesA[2,1] <-CuantilesData[7]; CuantilesA[2,2] <-CuantilesData[10]
CuantilesA[3,1] <-CuantilesData[6]; CuantilesA[3,2] <-CuantilesData[11]
CuantilesA[4,1] <-CuantilesData[5]; CuantilesA[4,2] <-CuantilesData[12]
CuantilesA[5,1] <-CuantilesData[4]; CuantilesA[5,2] <-CuantilesData[13]
CuantilesA[6,1] <-CuantilesData[3]; CuantilesA[6,2] <-CuantilesData[14]
CuantilesA[7,1] <-CuantilesData[2]; CuantilesA[7,2] <-CuantilesData[15]
CuantilesA[8,1] <-CuantilesData[1]; CuantilesA[8,2] <-CuantilesData[16]
CuantilesD[1,1] <-Cuantilillos[8];  CuantilesD[1,2] <-Cuantilillos[9]
CuantilesD[2,1] <-Cuantilillos[7];  CuantilesD[2,2] <-Cuantilillos[10]
CuantilesD[3,1] <-Cuantilillos[6];  CuantilesD[3,2] <-Cuantilillos[11]
CuantilesD[4,1] <-Cuantilillos[5];  CuantilesD[4,2] <-Cuantilillos[12]
CuantilesD[5,1] <-Cuantilillos[4];  CuantilesD[5,2] <-Cuantilillos[13]
CuantilesD[6,1] <-Cuantilillos[3];  CuantilesD[6,2] <-Cuantilillos[14]
CuantilesD[7,1] <-Cuantilillos[2];  CuantilesD[7,2] <-Cuantilillos[15]
CuantilesD[8,1] <-Cuantilillos[1];  CuantilesD[8,2] <-Cuantilillos[16]

print(CuantilesD)
##        LimInf    LimSup
## 65 -0.9345893 0.9345893
## 70 -1.0364334 1.0364334
## 75 -1.1503494 1.1503494
## 80 -1.2815516 1.2815516
## 85 -1.4395315 1.4395315
## 90 -1.6448536 1.6448536
## 95 -1.9599640 1.9599640
## 99 -2.5758293 2.5758293
print(CuantilesA)
##        LimInf    LimSup
## 65 -0.8936351 0.8467827
## 70 -0.8936351 0.9884290
## 75 -1.0112465 1.1344069
## 80 -1.0836965 1.3626690
## 85 -1.1697595 1.6777371
## 90 -1.2757625 2.1128235
## 95 -1.6121684 2.5115181
## 99 -1.9690098 3.0656774
col_sequence <- rainbow(n = 7, alpha = 0.35, start = 0, end = 1)
par(mfrow=c(2,1))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2 pEhEx - DATA (sample 2)', lty = 9)

abline(v=CuantilesA[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesA[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesA[2,1], lty=2, col="darkblue"); # 70% INFERIOR
abline(v=CuantilesA[2,2], lty=2, col="darkblue") # 70% SUPERIOR
abline(v=CuantilesA[3,1], lty=2, col="aquamarine4"); # 75% INFERIOR
abline(v=CuantilesA[3,2], lty=2, col="aquamarine4"); # 75% SUPERIOR
abline(v=CuantilesA[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesA[4,2], lty=2, col="green"); # 80% SUPERIOR
abline(v=CuantilesA[5,1], lty=2, col="brown"); # 85% INFERIOR
abline(v=CuantilesA[5,2], lty=2, col="brown"); # 85% SUPERIOR
abline(v=CuantilesA[6,1], lty=2, col="red");  # 90% INFERIOR
abline(v=CuantilesA[6,2], lty=2, col="red"); # 90% SUPERIOR
abline(v=CuantilesA[7,1], lty=2, col="blue");  # 95% INFERIOR
abline(v=CuantilesA[7,2], lty=2, col="blue"); # 95% SUPERIOR
abline(v=CuantilesA[8,1], lty=2, col="orange");  # 99% INFERIOR
abline(v=CuantilesA[8,2], lty=2, col="orange"); # 99% SUPERIOR
legend("topright",
       legend=c("65%","70%","75%","80%","85%","90%","95%","99%"), 
       pch=c(1,2,3,4,5,6,7,8),
       col=c("darkgoldenrod4","darkblue","aquamarine4",
             "green", "brown","red","blue","orange"))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2   pEhEx - ADJUSTED (sample 2)', lty=9)

abline(v=CuantilesD[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesD[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesD[2,1], lty=2, col="darkblue");  # 70% INFERIOR
abline(v=CuantilesD[2,2], lty=2, col="darkblue"); # 70% SUPERIOR
abline(v=CuantilesD[3,1], lty=2, col="aquamarine4");  # 75% INFERIOR
abline(v=CuantilesD[3,2], lty=2, col="aquamarine4"); # 75% SUPERIOR
abline(v=CuantilesD[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesD[4,2], lty=2, col="green"); # 80% SUPERIOR
abline(v=CuantilesD[5,1], lty=2, col="brown");  # 85% INFERIOR
abline(v=CuantilesD[5,2], lty=2, col="brown"); # 85% SUPERIOR
abline(v=CuantilesD[6,1], lty=2, col="red");  # 90% INFERIOR
abline(v=CuantilesD[6,2], lty=2, col="red"); # 90% SUPERIOR
abline(v=CuantilesD[7,1], lty=2, col="blue");  # 95% INFERIOR
abline(v=CuantilesD[7,2], lty=2, col="blue"); # 95% SUPERIOR
abline(v=CuantilesD[8,1], lty=2, col="orange");  # 99% INFERIOR
abline(v=CuantilesD[8,2], lty=2, col="orange"); # 99% SUPERIOR
legend("topright",
       legend=c("65%","70%","75%","80%","85%","90%","95%","99%"), 
       pch=c(1,2,3,4,5,6,7,8),
       col=c("darkgoldenrod4","darkblue","aquamarine4",
             "green", "brown","red","blue","orange"))

Grafica Cuantiles del \(65\%\) y \(80\%\)

par(mfrow=c(2,1))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2  pEhEx - DATA (sample 2)', lty=9)

abline(v=CuantilesA[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesA[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesA[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesA[4,2], lty=2, col="green"); # 80% SUPERIOR
legend("topright",legend=c("65%","80%"), 
       pch=c(1,2),col=c("darkgoldenrod4","green"))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2  BaseMean - ADJUSTED (sample 2)', lty=9)

abline(v=CuantilesD[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesD[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesD[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesD[4,2], lty=2, col="green"); # 80% SUPERIOR
legend("topright",legend=c("65%","80%"), 
       pch=c(1,2),col=c("darkgoldenrod4","green"))

Grafica Cuantiles del \(70\%\) y \(85\%\)

par(mfrow=c(2,1))
hist(tst, breaks = nbreaks, col = col_sequence, 
     main = 'Normalized Log2 pEhEx - DATA (sample 2)', lty=9)

abline(v=CuantilesA[2,1], lty=2, col="darkblue"); 
abline(v=CuantilesA[2,2], lty=2, col="darkblue")
abline(v=CuantilesA[5,1], lty=2, col="brown"); 
abline(v=CuantilesA[5,2], lty=2, col="brown")
legend("topright",legend=c("70%","85%"),
       pch=c(1,2),#3,4,5,6,7,8),
       col=c("brown"))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2  pEhEx - ADJUSTED (sample 2)', lty=9)

abline(v=CuantilesD[2,1], lty=2, col="darkblue"); 
abline(v=CuantilesD[2,2], lty=2, col="darkblue")
abline(v=CuantilesD[5,1], lty=2, col="brown"); 
abline(v=CuantilesD[5,2], lty=2, col="brown")
legend("topright",legend=c("70%","85%"),
       pch=c(1,2),#3,4,5,6,7,8),
       col=c("brown"))

Muestra 3

log2sample3 <- data5$log2sample3; head(mean(log2sample3)); head(sd(log2sample3))
## [1] 6.152478
## [1] 3.158081
head(log2sample3,5)
## [1]  8.904774  5.314000 11.052005  9.365349  5.680082
ndata5    <- length(log2sample3)
hist(log2sample3, breaks = nbreaks, col= rainbow(25,0.3), 
     main = 'Log2 sample2')

meanlog2sample3 <- mean(log2sample3); head(meanlog2sample3)
## [1] 6.152478
StdDevlog2sample3 <- sd(log2sample3); head(StdDevlog2sample3)
## [1] 3.158081
Normlog2sample3 <- (log2sample3-meanlog2sample3)/StdDevlog2sample3; head(Normlog2sample3)
## [1]  0.87150918 -0.26550237  1.55142547  1.01734948 -0.14958317 -0.05720518
tst<- Normlog2sample3
hist(tst, breaks = nbreaks, col= 1:5, 
     main = 'Normalized Log2sample1',
     xlab='pEhEx1',
     ylab= 'Frequency pEhEx')

Ajustando Modelos

fw1<-fitdist(tst, "norm")
plotdist(tst, histo = TRUE, demp = TRUE)

nnorm.f <- fitdist(tst,"norm")
summary(nnorm.f)
## Fitting of the distribution ' norm ' by maximum likelihood 
## Parameters : 
##          estimate Std. Error
## mean 1.292160e-16 0.01425665
## sd   9.998983e-01 0.01008093
## Loglikelihood:  -6979.259   AIC:  13962.52   BIC:  13975.52 
## Correlation matrix:
##              mean           sd
## mean  1.00000e+00 -3.26782e-11
## sd   -3.26782e-11  1.00000e+00
par(mfrow=c(2,2))
denscomp(nnorm.f,legendtext = 'Dist Normal')
qqcomp(nnorm.f,legendtext = 'Dist Normal')
cdfcomp(nnorm.f,legendtext = 'Dist Normal')
ppcomp(nnorm.f,legendtext = 'Dist Normal')

probs <- c();
probs[8] = 0.175;  probs[9] = 0.825; 
probs[7] = 0.15;   probs[10] = 0.85;   
probs[6] = 0.125;  probs[11] = 0.875; 
probs[5] = 0.1;    probs[12] = 0.9;    
probs[4] = 0.075;  probs[13] = 0.925; 
probs[3] = 0.05;   probs[14] = 0.95;   
probs[2] = 0.025;  probs[15] = 0.975; 
probs[1] = 0.005;  probs[16] = 0.995;  
CuantilesData <- quantile(tst,prob = probs)
CuantilesModel <- qnorm(probs, mean=0, sd=1)
Cuantilillos <- t(CuantilesModel)
colnames(Cuantilillos) <- c('0.5%','2.5%','5%','7.5%',
                            '10%','12.5%','15%','17.5%',
                            '82.5%','85%','87.5%','90%',
                            '92.5%','95%','97.5%','99.5%')
Cuantilillos <- t(Cuantilillos)
colnames(Cuantilillos) <- c('Cuantiles Ajuste')
print(Cuantilillos)
##       Cuantiles Ajuste
## 0.5%        -2.5758293
## 2.5%        -1.9599640
## 5%          -1.6448536
## 7.5%        -1.4395315
## 10%         -1.2815516
## 12.5%       -1.1503494
## 15%         -1.0364334
## 17.5%       -0.9345893
## 82.5%        0.9345893
## 85%          1.0364334
## 87.5%        1.1503494
## 90%          1.2815516
## 92.5%        1.4395315
## 95%          1.6448536
## 97.5%        1.9599640
## 99.5%        2.5758293
CuantilesA <- matrix(0,8,2)
CuantilesD <- matrix(0,8,2)
colnames(CuantilesA) <- c('LimInf','LimSup')
colnames(CuantilesD) <- c('LimInf','LimSup')
rownames(CuantilesA) <- c('65','70','75','80','85','90','95','99')
rownames(CuantilesD) <- c('65','70','75','80','85','90','95','99')
CuantilesA[1,1] <-CuantilesData[8]; CuantilesA[1,2] <-CuantilesData[9]
CuantilesA[2,1] <-CuantilesData[7]; CuantilesA[2,2] <-CuantilesData[10]
CuantilesA[3,1] <-CuantilesData[6]; CuantilesA[3,2] <-CuantilesData[11]
CuantilesA[4,1] <-CuantilesData[5]; CuantilesA[4,2] <-CuantilesData[12]
CuantilesA[5,1] <-CuantilesData[4]; CuantilesA[5,2] <-CuantilesData[13]
CuantilesA[6,1] <-CuantilesData[3]; CuantilesA[6,2] <-CuantilesData[14]
CuantilesA[7,1] <-CuantilesData[2]; CuantilesA[7,2] <-CuantilesData[15]
CuantilesA[8,1] <-CuantilesData[1]; CuantilesA[8,2] <-CuantilesData[16]
CuantilesD[1,1] <-Cuantilillos[8];  CuantilesD[1,2] <-Cuantilillos[9]
CuantilesD[2,1] <-Cuantilillos[7];  CuantilesD[2,2] <-Cuantilillos[10]
CuantilesD[3,1] <-Cuantilillos[6];  CuantilesD[3,2] <-Cuantilillos[11]
CuantilesD[4,1] <-Cuantilillos[5];  CuantilesD[4,2] <-Cuantilillos[12]
CuantilesD[5,1] <-Cuantilillos[4];  CuantilesD[5,2] <-Cuantilillos[13]
CuantilesD[6,1] <-Cuantilillos[3];  CuantilesD[6,2] <-Cuantilillos[14]
CuantilesD[7,1] <-Cuantilillos[2];  CuantilesD[7,2] <-Cuantilillos[15]
CuantilesD[8,1] <-Cuantilillos[1];  CuantilesD[8,2] <-Cuantilillos[16]
print(CuantilesD)
##        LimInf    LimSup
## 65 -0.9345893 0.9345893
## 70 -1.0364334 1.0364334
## 75 -1.1503494 1.1503494
## 80 -1.2815516 1.2815516
## 85 -1.4395315 1.4395315
## 90 -1.6448536 1.6448536
## 95 -1.9599640 1.9599640
## 99 -2.5758293 2.5758293
print(CuantilesA)
##        LimInf    LimSup
## 65 -0.8811857 0.8359346
## 70 -0.9495345 0.9646852
## 75 -1.0299396 1.1073525
## 80 -1.1855441 1.3225658
## 85 -1.2519401 1.6081997
## 90 -1.4233556 2.0263238
## 95 -1.7009100 2.4893700
## 99 -1.9481699 3.1039165

CreaciĂ³n de histogramas

col_sequence <- rainbow(n = 7, alpha = 0.35, start = 0, end = 1)
par(mfrow=c(2,1))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2 pEhEx - DATA (sample 3)', lty = 9)

abline(v=CuantilesA[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesA[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesA[2,1], lty=2, col="darkblue"); # 70% INFERIOR
abline(v=CuantilesA[2,2], lty=2, col="darkblue") # 70% SUPERIOR
abline(v=CuantilesA[3,1], lty=2, col="aquamarine4"); # 75% INFERIOR
abline(v=CuantilesA[3,2], lty=2, col="aquamarine4"); # 75% SUPERIOR
abline(v=CuantilesA[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesA[4,2], lty=2, col="green"); # 80% SUPERIOR
abline(v=CuantilesA[5,1], lty=2, col="brown"); # 85% INFERIOR
abline(v=CuantilesA[5,2], lty=2, col="brown"); # 85% SUPERIOR
abline(v=CuantilesA[6,1], lty=2, col="red");  # 90% INFERIOR
abline(v=CuantilesA[6,2], lty=2, col="red"); # 90% SUPERIOR
abline(v=CuantilesA[7,1], lty=2, col="blue");  # 95% INFERIOR
abline(v=CuantilesA[7,2], lty=2, col="blue"); # 95% SUPERIOR
abline(v=CuantilesA[8,1], lty=2, col="orange");  # 99% INFERIOR
abline(v=CuantilesA[8,2], lty=2, col="orange"); # 99% SUPERIOR
legend("topright",
       legend=c("65%","70%","75%","80%","85%","90%","95%","99%"), 
       pch=c(1,2,3,4,5,6,7,8),
       col=c("darkgoldenrod4","darkblue","aquamarine4",
             "green", "brown","red","blue","orange"))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2   pEhEx - ADJUSTED (sample 3)', lty=9)

abline(v=CuantilesD[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesD[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesD[2,1], lty=2, col="darkblue");  # 70% INFERIOR
abline(v=CuantilesD[2,2], lty=2, col="darkblue"); # 70% SUPERIOR
abline(v=CuantilesD[3,1], lty=2, col="aquamarine4");  # 75% INFERIOR
abline(v=CuantilesD[3,2], lty=2, col="aquamarine4"); # 75% SUPERIOR
abline(v=CuantilesD[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesD[4,2], lty=2, col="green"); # 80% SUPERIOR
abline(v=CuantilesD[5,1], lty=2, col="brown");  # 85% INFERIOR
abline(v=CuantilesD[5,2], lty=2, col="brown"); # 85% SUPERIOR
abline(v=CuantilesD[6,1], lty=2, col="red");  # 90% INFERIOR
abline(v=CuantilesD[6,2], lty=2, col="red"); # 90% SUPERIOR
abline(v=CuantilesD[7,1], lty=2, col="blue");  # 95% INFERIOR
abline(v=CuantilesD[7,2], lty=2, col="blue"); # 95% SUPERIOR
abline(v=CuantilesD[8,1], lty=2, col="orange");  # 99% INFERIOR
abline(v=CuantilesD[8,2], lty=2, col="orange"); # 99% SUPERIOR
legend("topright",
       legend=c("65%","70%","75%","80%","85%","90%","95%","99%"), 
       pch=c(1,2,3,4,5,6,7,8),
       col=c("darkgoldenrod4","darkblue","aquamarine4",
             "green", "brown","red","blue","orange"))

Grafica Cuantiles del \(65\%\) y \(80\%\)

par(mfrow=c(2,1))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2  pEhEx - DATA (sample 2)', lty=9)

abline(v=CuantilesA[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesA[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesA[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesA[4,2], lty=2, col="green"); # 80% SUPERIOR
legend("topright",legend=c("65%","80%"), 
       pch=c(1,2),col=c("darkgoldenrod4","green"))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2  BaseMean - ADJUSTED (sample 3)', lty=9)

abline(v=CuantilesD[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesD[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesD[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesD[4,2], lty=2, col="green"); # 80% SUPERIOR
legend("topright",legend=c("65%","80%"), 
       pch=c(1,2),col=c("darkgoldenrod4","green"))

Grafica Cuantiles del \(70\%\) y \(85\%\)

par(mfrow=c(2,1))
hist(tst, breaks = nbreaks, col = col_sequence, 
     main = 'Normalized Log2 pEhEx - DATA (sample 3)', lty=9)

abline(v=CuantilesA[2,1], lty=2, col="darkblue"); 
abline(v=CuantilesA[2,2], lty=2, col="darkblue")
abline(v=CuantilesA[5,1], lty=2, col="brown"); 
abline(v=CuantilesA[5,2], lty=2, col="brown")
legend("topright",legend=c("70%","85%"),
       pch=c(1,2),#3,4,5,6,7,8),
       col=c("brown"))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2  pEhEx - ADJUSTED (sample 3)', lty=9)

abline(v=CuantilesD[2,1], lty=2, col="darkblue"); 
abline(v=CuantilesD[2,2], lty=2, col="darkblue")
abline(v=CuantilesD[5,1], lty=2, col="brown"); 
abline(v=CuantilesD[5,2], lty=2, col="brown")
legend("topright",legend=c("70%","85%"),
       pch=c(1,2),#3,4,5,6,7,8),
       col=c("brown"))

Muestra pEhExvsumasM1

log2vsumasM1 <- data5$log2samplevsumasM1; head(mean(log2vsumasM1)); head(sd(log2vsumasM1))
## [1] 6.221901
## [1] 2.955459
head(log2vsumasM1,5)
## [1] 5.749419 6.889844 9.540979 9.462180 6.904367
ndata5    <- length(log2vsumasM1)
hist(log2vsumasM1, breaks = nbreaks, col= rainbow(25,0.3), 
     main = 'Log2vsumasM1')

meanlog2vsumasM1 <- mean(log2vsumasM1); head(meanlog2vsumasM1)
## [1] 6.221901
StdDevlog2vsumasM1 <- sd(log2vsumasM1); head(StdDevlog2vsumasM1)
## [1] 2.955459
Normlog2vsumasM1 <- (log2vsumasM1-meanlog2vsumasM1)/StdDevlog2vsumasM1; head(Normlog2vsumasM1)
## [1] -0.1598676  0.2260030  1.1230332  1.0963710  0.2309173 -0.1181501
tst<- Normlog2vsumasM1

Primer histograma

hist(tst, breaks = nbreaks, col= 1:5, 
     main = 'Normalized Log2vsumasM1',
     xlab='pEhEx1',
     ylab= 'Frequency pEhEx')

Ajustando modelo

fw1<-fitdist(tst, "norm")
plotdist(tst, histo = TRUE, demp = TRUE)

nnorm.f <- fitdist(tst,"norm")
summary(nnorm.f)
## Fitting of the distribution ' norm ' by maximum likelihood 
## Parameters : 
##           estimate Std. Error
## mean -9.853035e-17 0.01425665
## sd    9.998983e-01 0.01008093
## Loglikelihood:  -6979.259   AIC:  13962.52   BIC:  13975.52 
## Correlation matrix:
##      mean sd
## mean    1  0
## sd      0  1
par(mfrow=c(2,2))
denscomp(nnorm.f,legendtext = 'Dist Normal')
qqcomp(nnorm.f,legendtext = 'Dist Normal')
cdfcomp(nnorm.f,legendtext = 'Dist Normal')
ppcomp(nnorm.f,legendtext = 'Dist Normal')

probs <- c();
probs[8] = 0.175;  probs[9] = 0.825; 
probs[7] = 0.15;   probs[10] = 0.85;   
probs[6] = 0.125;  probs[11] = 0.875; 
probs[5] = 0.1;    probs[12] = 0.9;    
probs[4] = 0.075;  probs[13] = 0.925; 
probs[3] = 0.05;   probs[14] = 0.95;   
probs[2] = 0.025;  probs[15] = 0.975; 
probs[1] = 0.005;  probs[16] = 0.995;  
CuantilesData <- quantile(tst,prob = probs)
CuantilesModel <- qnorm(probs, mean=0, sd=1)
Cuantilillos <- t(CuantilesModel)
colnames(Cuantilillos) <- c('0.5%','2.5%','5%','7.5%',
                            '10%','12.5%','15%','17.5%',
                            '82.5%','85%','87.5%','90%',
                            '92.5%','95%','97.5%','99.5%')
Cuantilillos <- t(Cuantilillos)
colnames(Cuantilillos) <- c('Cuantiles Ajuste')
print(Cuantilillos)
##       Cuantiles Ajuste
## 0.5%        -2.5758293
## 2.5%        -1.9599640
## 5%          -1.6448536
## 7.5%        -1.4395315
## 10%         -1.2815516
## 12.5%       -1.1503494
## 15%         -1.0364334
## 17.5%       -0.9345893
## 82.5%        0.9345893
## 85%          1.0364334
## 87.5%        1.1503494
## 90%          1.2815516
## 92.5%        1.4395315
## 95%          1.6448536
## 97.5%        1.9599640
## 99.5%        2.5758293
CuantilesA <- matrix(0,8,2)
CuantilesD <- matrix(0,8,2)
colnames(CuantilesA) <- c('LimInf','LimSup')
colnames(CuantilesD) <- c('LimInf','LimSup')
rownames(CuantilesA) <- c('65','70','75','80','85','90','95','99')
rownames(CuantilesD) <- c('65','70','75','80','85','90','95','99')
CuantilesA[1,1] <-CuantilesData[8]; CuantilesA[1,2] <-CuantilesData[9]
CuantilesA[2,1] <-CuantilesData[7]; CuantilesA[2,2] <-CuantilesData[10]
CuantilesA[3,1] <-CuantilesData[6]; CuantilesA[3,2] <-CuantilesData[11]
CuantilesA[4,1] <-CuantilesData[5]; CuantilesA[4,2] <-CuantilesData[12]
CuantilesA[5,1] <-CuantilesData[4]; CuantilesA[5,2] <-CuantilesData[13]
CuantilesA[6,1] <-CuantilesData[3]; CuantilesA[6,2] <-CuantilesData[14]
CuantilesA[7,1] <-CuantilesData[2]; CuantilesA[7,2] <-CuantilesData[15]
CuantilesA[8,1] <-CuantilesData[1]; CuantilesA[8,2] <-CuantilesData[16]
CuantilesD[1,1] <-Cuantilillos[8];  CuantilesD[1,2] <-Cuantilillos[9]
CuantilesD[2,1] <-Cuantilillos[7];  CuantilesD[2,2] <-Cuantilillos[10]
CuantilesD[3,1] <-Cuantilillos[6];  CuantilesD[3,2] <-Cuantilillos[11]
CuantilesD[4,1] <-Cuantilillos[5];  CuantilesD[4,2] <-Cuantilillos[12]
CuantilesD[5,1] <-Cuantilillos[4];  CuantilesD[5,2] <-Cuantilillos[13]
CuantilesD[6,1] <-Cuantilillos[3];  CuantilesD[6,2] <-Cuantilillos[14]
CuantilesD[7,1] <-Cuantilillos[2];  CuantilesD[7,2] <-Cuantilillos[15]
CuantilesD[8,1] <-Cuantilillos[1];  CuantilesD[8,2] <-Cuantilillos[16]
print(CuantilesD)
##        LimInf    LimSup
## 65 -0.9345893 0.9345893
## 70 -1.0364334 1.0364334
## 75 -1.1503494 1.1503494
## 80 -1.2815516 1.2815516
## 85 -1.4395315 1.4395315
## 90 -1.6448536 1.6448536
## 95 -1.9599640 1.9599640
## 99 -2.5758293 2.5758293
print(CuantilesA)
##        LimInf    LimSup
## 65 -0.8531906 0.8236011
## 70 -0.9004672 0.9332295
## 75 -0.9528182 1.1118993
## 80 -0.9528182 1.3280269
## 85 -1.0781339 1.6545578
## 90 -1.1553694 2.1045904
## 95 -1.3603164 2.6974724
## 99 -2.1052233 3.3025359

Histogramas

col_sequence <- rainbow(n = 7, alpha = 0.35, start = 0, end = 1)
par(mfrow=c(2,1))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2vsumasM1 pEhEx - DATA', lty = 9)

abline(v=CuantilesA[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesA[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesA[2,1], lty=2, col="darkblue"); # 70% INFERIOR
abline(v=CuantilesA[2,2], lty=2, col="darkblue") # 70% SUPERIOR
abline(v=CuantilesA[3,1], lty=2, col="aquamarine4"); # 75% INFERIOR
abline(v=CuantilesA[3,2], lty=2, col="aquamarine4"); # 75% SUPERIOR
abline(v=CuantilesA[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesA[4,2], lty=2, col="green"); # 80% SUPERIOR
abline(v=CuantilesA[5,1], lty=2, col="brown"); # 85% INFERIOR
abline(v=CuantilesA[5,2], lty=2, col="brown"); # 85% SUPERIOR
abline(v=CuantilesA[6,1], lty=2, col="red");  # 90% INFERIOR
abline(v=CuantilesA[6,2], lty=2, col="red"); # 90% SUPERIOR
abline(v=CuantilesA[7,1], lty=2, col="blue");  # 95% INFERIOR
abline(v=CuantilesA[7,2], lty=2, col="blue"); # 95% SUPERIOR
abline(v=CuantilesA[8,1], lty=2, col="orange");  # 99% INFERIOR
abline(v=CuantilesA[8,2], lty=2, col="orange"); # 99% SUPERIOR
legend("topright",
       legend=c("65%","70%","75%","80%","85%","90%","95%","99%"), 
       pch=c(1,2,3,4,5,6,7,8),
       col=c("darkgoldenrod4","darkblue","aquamarine4",
             "green", "brown","red","blue","orange"))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2vsumasM1   pEhEx - ADJUSTED', lty=9)

abline(v=CuantilesD[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesD[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesD[2,1], lty=2, col="darkblue");  # 70% INFERIOR
abline(v=CuantilesD[2,2], lty=2, col="darkblue"); # 70% SUPERIOR
abline(v=CuantilesD[3,1], lty=2, col="aquamarine4");  # 75% INFERIOR
abline(v=CuantilesD[3,2], lty=2, col="aquamarine4"); # 75% SUPERIOR
abline(v=CuantilesD[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesD[4,2], lty=2, col="green"); # 80% SUPERIOR
abline(v=CuantilesD[5,1], lty=2, col="brown");  # 85% INFERIOR
abline(v=CuantilesD[5,2], lty=2, col="brown"); # 85% SUPERIOR
abline(v=CuantilesD[6,1], lty=2, col="red");  # 90% INFERIOR
abline(v=CuantilesD[6,2], lty=2, col="red"); # 90% SUPERIOR
abline(v=CuantilesD[7,1], lty=2, col="blue");  # 95% INFERIOR
abline(v=CuantilesD[7,2], lty=2, col="blue"); # 95% SUPERIOR
abline(v=CuantilesD[8,1], lty=2, col="orange");  # 99% INFERIOR
abline(v=CuantilesD[8,2], lty=2, col="orange"); # 99% SUPERIOR
legend("topright",
       legend=c("65%","70%","75%","80%","85%","90%","95%","99%"), 
       pch=c(1,2,3,4,5,6,7,8),
       col=c("darkgoldenrod4","darkblue","aquamarine4",
             "green", "brown","red","blue","orange"))

Grafica Cuantiles del \(65\%\) y \(80\%\)

par(mfrow=c(2,1))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2vsumasM1  pEhEx - DATA', lty=9)

abline(v=CuantilesA[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesA[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesA[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesA[4,2], lty=2, col="green"); # 80% SUPERIOR
legend("topright",legend=c("65%","80%"), 
       pch=c(1,2),col=c("darkgoldenrod4","green"))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2vsumasM1  - ADJUSTED', lty=9)

abline(v=CuantilesD[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesD[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesD[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesD[4,2], lty=2, col="green"); # 80% SUPERIOR
legend("topright",legend=c("65%","80%"), 
       pch=c(1,2),col=c("darkgoldenrod4","green"))

Grafica Cuantiles del \(70\%\) y \(85\%\)

par(mfrow=c(2,1))
hist(tst, breaks = nbreaks, col = col_sequence, 
     main = 'Normalized Log2vsumasM1 pEhEx - DATA', lty=9)

abline(v=CuantilesA[2,1], lty=2, col="darkblue"); 
abline(v=CuantilesA[2,2], lty=2, col="darkblue")
abline(v=CuantilesA[5,1], lty=2, col="brown"); 
abline(v=CuantilesA[5,2], lty=2, col="brown")
legend("topright",legend=c("70%","85%"),
       pch=c(1,2),#3,4,5,6,7,8),
       col=c("brown"))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2vsumasM1  pEhEx - ADJUSTED', lty=9)

abline(v=CuantilesD[2,1], lty=2, col="darkblue"); 
abline(v=CuantilesD[2,2], lty=2, col="darkblue")
abline(v=CuantilesD[5,1], lty=2, col="brown"); 
abline(v=CuantilesD[5,2], lty=2, col="brown")
legend("topright",legend=c("70%","85%"),
       pch=c(1,2),#3,4,5,6,7,8),
       col=c("brown"))

Muestra pEhExvsumasM2

log2vsumasM2 <- data5$log2samplevsumasM2; head(mean(log2vsumasM2)); head(sd(log2vsumasM2))
## [1] 6.255602
## [1] 2.72471
head(log2vsumasM2,5)
## [1]  8.745886  8.538024 10.187333  7.740125  6.603402

Primer Histograma

ndata5    <- length(log2vsumasM2)
hist(log2vsumasM2, breaks = nbreaks, col= rainbow(25,0.3), 
     main = 'Log2vsumasM2')

meanlog2vsumasM2 <- mean(log2vsumasM2); head(meanlog2vsumasM2)
## [1] 6.255602
StdDevlog2vsumasM2 <- sd(log2vsumasM2); head(StdDevlog2vsumasM2)
## [1] 2.72471
Normlog2vsumasM2 <- (log2vsumasM2-meanlog2vsumasM2)/StdDevlog2vsumasM2; head(Normlog2vsumasM2)
## [1]  0.9139630  0.8376754  1.4429905  0.5448372  0.1276466 -0.1223761
tst<- Normlog2vsumasM2

** Segundo Histograma**

hist(tst, breaks = nbreaks, col= 1:5, 
     main = 'Normalized Log2vsumasM2',
     xlab='pEhEx1',
     ylab= 'Frequency pEhEx')

Ajustando modelo

fw1<-fitdist(tst, "norm")
plotdist(tst, histo = TRUE, demp = TRUE)

nnorm.f <- fitdist(tst,"norm")
summary(nnorm.f)
## Fitting of the distribution ' norm ' by maximum likelihood 
## Parameters : 
##          estimate Std. Error
## mean 1.626526e-16 0.01425665
## sd   9.998983e-01 0.01008093
## Loglikelihood:  -6979.259   AIC:  13962.52   BIC:  13975.52 
## Correlation matrix:
##      mean sd
## mean    1  0
## sd      0  1
par(mfrow=c(2,2))
denscomp(nnorm.f,legendtext = 'Dist Normal')
qqcomp(nnorm.f,legendtext = 'Dist Normal')
cdfcomp(nnorm.f,legendtext = 'Dist Normal')
ppcomp(nnorm.f,legendtext = 'Dist Normal')

CĂ¡lculo de cuantiles

probs <- c();
probs[8] = 0.175;  probs[9] = 0.825; 
probs[7] = 0.15;   probs[10] = 0.85;   
probs[6] = 0.125;  probs[11] = 0.875; 
probs[5] = 0.1;    probs[12] = 0.9;    
probs[4] = 0.075;  probs[13] = 0.925; 
probs[3] = 0.05;   probs[14] = 0.95;   
probs[2] = 0.025;  probs[15] = 0.975; 
probs[1] = 0.005;  probs[16] = 0.995;  
CuantilesData <- quantile(tst,prob = probs)
CuantilesModel <- qnorm(probs, mean=0, sd=1)
Cuantilillos <- t(CuantilesModel)
colnames(Cuantilillos) <- c('0.5%','2.5%','5%','7.5%',
                            '10%','12.5%','15%','17.5%',
                            '82.5%','85%','87.5%','90%',
                            '92.5%','95%','97.5%','99.5%')
Cuantilillos <- t(Cuantilillos)
colnames(Cuantilillos) <- c('Cuantiles Ajuste')
print(Cuantilillos)
##       Cuantiles Ajuste
## 0.5%        -2.5758293
## 2.5%        -1.9599640
## 5%          -1.6448536
## 7.5%        -1.4395315
## 10%         -1.2815516
## 12.5%       -1.1503494
## 15%         -1.0364334
## 17.5%       -0.9345893
## 82.5%        0.9345893
## 85%          1.0364334
## 87.5%        1.1503494
## 90%          1.2815516
## 92.5%        1.4395315
## 95%          1.6448536
## 97.5%        1.9599640
## 99.5%        2.5758293
CuantilesA <- matrix(0,8,2)
CuantilesD <- matrix(0,8,2)
colnames(CuantilesA) <- c('LimInf','LimSup')
colnames(CuantilesD) <- c('LimInf','LimSup')
rownames(CuantilesA) <- c('65','70','75','80','85','90','95','99')
rownames(CuantilesD) <- c('65','70','75','80','85','90','95','99')
CuantilesA[1,1] <-CuantilesData[8]; CuantilesA[1,2] <-CuantilesData[9]
CuantilesA[2,1] <-CuantilesData[7]; CuantilesA[2,2] <-CuantilesData[10]
CuantilesA[3,1] <-CuantilesData[6]; CuantilesA[3,2] <-CuantilesData[11]
CuantilesA[4,1] <-CuantilesData[5]; CuantilesA[4,2] <-CuantilesData[12]
CuantilesA[5,1] <-CuantilesData[4]; CuantilesA[5,2] <-CuantilesData[13]
CuantilesA[6,1] <-CuantilesData[3]; CuantilesA[6,2] <-CuantilesData[14]
CuantilesA[7,1] <-CuantilesData[2]; CuantilesA[7,2] <-CuantilesData[15]
CuantilesA[8,1] <-CuantilesData[1]; CuantilesA[8,2] <-CuantilesData[16]
CuantilesD[1,1] <-Cuantilillos[8];  CuantilesD[1,2] <-Cuantilillos[9]
CuantilesD[2,1] <-Cuantilillos[7];  CuantilesD[2,2] <-Cuantilillos[10]
CuantilesD[3,1] <-Cuantilillos[6];  CuantilesD[3,2] <-Cuantilillos[11]
CuantilesD[4,1] <-Cuantilillos[5];  CuantilesD[4,2] <-Cuantilillos[12]
CuantilesD[5,1] <-Cuantilillos[4];  CuantilesD[5,2] <-Cuantilillos[13]
CuantilesD[6,1] <-Cuantilillos[3];  CuantilesD[6,2] <-Cuantilillos[14]
CuantilesD[7,1] <-Cuantilillos[2];  CuantilesD[7,2] <-Cuantilillos[15]
CuantilesD[8,1] <-Cuantilillos[1];  CuantilesD[8,2] <-Cuantilillos[16]
print(CuantilesD)
##        LimInf    LimSup
## 65 -0.9345893 0.9345893
## 70 -1.0364334 1.0364334
## 75 -1.1503494 1.1503494
## 80 -1.2815516 1.2815516
## 85 -1.4395315 1.4395315
## 90 -1.6448536 1.6448536
## 95 -1.9599640 1.9599640
## 99 -2.5758293 2.5758293
print(CuantilesA)
##        LimInf    LimSup
## 65 -0.8779788 0.8200071
## 70 -0.9297388 0.9699103
## 75 -0.9871123 1.1247423
## 80 -1.0514667 1.3424300
## 85 -1.1247392 1.7188713
## 90 -1.2098054 2.1841922
## 95 -1.3112040 2.5940822
## 99 -1.8422575 3.1504680

Histogramas

col_sequence <- rainbow(n = 7, alpha = 0.35, start = 0, end = 1)
par(mfrow=c(2,1))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2vsumasM2 pEhEx - DATA', lty = 9)

abline(v=CuantilesA[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesA[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesA[2,1], lty=2, col="darkblue"); # 70% INFERIOR
abline(v=CuantilesA[2,2], lty=2, col="darkblue") # 70% SUPERIOR
abline(v=CuantilesA[3,1], lty=2, col="aquamarine4"); # 75% INFERIOR
abline(v=CuantilesA[3,2], lty=2, col="aquamarine4"); # 75% SUPERIOR
abline(v=CuantilesA[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesA[4,2], lty=2, col="green"); # 80% SUPERIOR
abline(v=CuantilesA[5,1], lty=2, col="brown"); # 85% INFERIOR
abline(v=CuantilesA[5,2], lty=2, col="brown"); # 85% SUPERIOR
abline(v=CuantilesA[6,1], lty=2, col="red");  # 90% INFERIOR
abline(v=CuantilesA[6,2], lty=2, col="red"); # 90% SUPERIOR
abline(v=CuantilesA[7,1], lty=2, col="blue");  # 95% INFERIOR
abline(v=CuantilesA[7,2], lty=2, col="blue"); # 95% SUPERIOR
abline(v=CuantilesA[8,1], lty=2, col="orange");  # 99% INFERIOR
abline(v=CuantilesA[8,2], lty=2, col="orange"); # 99% SUPERIOR
legend("topright",
       legend=c("65%","70%","75%","80%","85%","90%","95%","99%"), 
       pch=c(1,2,3,4,5,6,7,8),
       col=c("darkgoldenrod4","darkblue","aquamarine4",
             "green", "brown","red","blue","orange"))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2vsumasM2   pEhEx - ADJUSTED', lty=9)

abline(v=CuantilesD[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesD[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesD[2,1], lty=2, col="darkblue");  # 70% INFERIOR
abline(v=CuantilesD[2,2], lty=2, col="darkblue"); # 70% SUPERIOR
abline(v=CuantilesD[3,1], lty=2, col="aquamarine4");  # 75% INFERIOR
abline(v=CuantilesD[3,2], lty=2, col="aquamarine4"); # 75% SUPERIOR
abline(v=CuantilesD[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesD[4,2], lty=2, col="green"); # 80% SUPERIOR
abline(v=CuantilesD[5,1], lty=2, col="brown");  # 85% INFERIOR
abline(v=CuantilesD[5,2], lty=2, col="brown"); # 85% SUPERIOR
abline(v=CuantilesD[6,1], lty=2, col="red");  # 90% INFERIOR
abline(v=CuantilesD[6,2], lty=2, col="red"); # 90% SUPERIOR
abline(v=CuantilesD[7,1], lty=2, col="blue");  # 95% INFERIOR
abline(v=CuantilesD[7,2], lty=2, col="blue"); # 95% SUPERIOR
abline(v=CuantilesD[8,1], lty=2, col="orange");  # 99% INFERIOR
abline(v=CuantilesD[8,2], lty=2, col="orange"); # 99% SUPERIOR
legend("topright",
       legend=c("65%","70%","75%","80%","85%","90%","95%","99%"), 
       pch=c(1,2,3,4,5,6,7,8),
       col=c("darkgoldenrod4","darkblue","aquamarine4",
             "green", "brown","red","blue","orange"))

Grafica Cuantiles del \(65\%\) y \(80\%\)

par(mfrow=c(2,1))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2vsumasM2  pEhEx - DATA', lty=9)

abline(v=CuantilesA[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesA[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesA[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesA[4,2], lty=2, col="green"); # 80% SUPERIOR
legend("topright",legend=c("65%","80%"), 
       pch=c(1,2),col=c("darkgoldenrod4","green"))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2vsumasM2  - ADJUSTED', lty=9)

abline(v=CuantilesD[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesD[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesD[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesD[4,2], lty=2, col="green"); # 80% SUPERIOR
legend("topright",legend=c("65%","80%"), 
       pch=c(1,2),col=c("darkgoldenrod4","green"))

Grafica Cuantiles del \(70\%\) y \(85\%\)

par(mfrow=c(2,1))
hist(tst, breaks = nbreaks, col = col_sequence, 
     main = 'Normalized Log2vsumasM2 pEhEx - DATA', lty=9)

abline(v=CuantilesA[2,1], lty=2, col="darkblue"); 
abline(v=CuantilesA[2,2], lty=2, col="darkblue")
abline(v=CuantilesA[5,1], lty=2, col="brown"); 
abline(v=CuantilesA[5,2], lty=2, col="brown")
legend("topright",legend=c("70%","85%"),
       pch=c(1,2),#3,4,5,6,7,8),
       col=c("brown"))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2vsumasM2  pEhEx - ADJUSTED', lty=9)

abline(v=CuantilesD[2,1], lty=2, col="darkblue"); 
abline(v=CuantilesD[2,2], lty=2, col="darkblue")
abline(v=CuantilesD[5,1], lty=2, col="brown"); 
abline(v=CuantilesD[5,2], lty=2, col="brown")
legend("topright",legend=c("70%","85%"),
       pch=c(1,2),#3,4,5,6,7,8),
       col=c("brown"))

Muestra pEhExvsumasM3

log2vsumasM3 <- data5$log2samplevsumasM3; head(mean(log2vsumasM3)); head(sd(log2vsumasM3))
## [1] 6.317092
## [1] 3.410882
head(log2vsumasM3,5)
## [1]  9.778619  5.482538 12.794720 10.625030  7.128629
ndata5    <- length(log2vsumasM3)

** Primer histograma**

hist(log2vsumasM3, breaks = nbreaks, col= rainbow(25,0.3), 
     main = 'Log2vsumasM3')

meanlog2vsumasM3 <- mean(log2vsumasM3); head(meanlog2vsumasM3)
## [1] 6.317092
StdDevlog2vsumasM3 <- sd(log2vsumasM3); head(StdDevlog2vsumasM3)
## [1] 3.410882
Normlog2vsumasM3 <- (log2vsumasM3-meanlog2vsumasM3)/StdDevlog2vsumasM3; head(Normlog2vsumasM3)
## [1]  1.0148481 -0.2446740  1.8991065  1.2629984  0.2379259 -0.3441057
tst<- Normlog2vsumasM3

Segundo histograma

hist(tst, breaks = nbreaks, col= 1:5, 
     main = 'Normalized Log2vsumasM3',
     xlab='pEhEx1',
     ylab= 'Frequency pEhEx')

Ajustando modelo

fw1<-fitdist(tst, "norm")
plotdist(tst, histo = TRUE, demp = TRUE)

nnorm.f <- fitdist(tst,"norm")
summary(nnorm.f)
## Fitting of the distribution ' norm ' by maximum likelihood 
## Parameters : 
##           estimate Std. Error
## mean -2.967667e-17 0.01425665
## sd    9.998983e-01 0.01008093
## Loglikelihood:  -6979.259   AIC:  13962.52   BIC:  13975.52 
## Correlation matrix:
##      mean sd
## mean    1  0
## sd      0  1
par(mfrow=c(2,2))
denscomp(nnorm.f,legendtext = 'Dist Normal')
qqcomp(nnorm.f,legendtext = 'Dist Normal')
cdfcomp(nnorm.f,legendtext = 'Dist Normal')
ppcomp(nnorm.f,legendtext = 'Dist Normal')

Calculo de cuantiles

probs <- c();
probs[8] = 0.175;  probs[9] = 0.825; 
probs[7] = 0.15;   probs[10] = 0.85;   
probs[6] = 0.125;  probs[11] = 0.875; 
probs[5] = 0.1;    probs[12] = 0.9;    
probs[4] = 0.075;  probs[13] = 0.925; 
probs[3] = 0.05;   probs[14] = 0.95;   
probs[2] = 0.025;  probs[15] = 0.975; 
probs[1] = 0.005;  probs[16] = 0.995;  
CuantilesData <- quantile(tst,prob = probs)
CuantilesModel <- qnorm(probs, mean=0, sd=1)
Cuantilillos <- t(CuantilesModel)
colnames(Cuantilillos) <- c('0.5%','2.5%','5%','7.5%',
                            '10%','12.5%','15%','17.5%',
                            '82.5%','85%','87.5%','90%',
                            '92.5%','95%','97.5%','99.5%')
Cuantilillos <- t(Cuantilillos)
colnames(Cuantilillos) <- c('Cuantiles Ajuste')
print(Cuantilillos)
##       Cuantiles Ajuste
## 0.5%        -2.5758293
## 2.5%        -1.9599640
## 5%          -1.6448536
## 7.5%        -1.4395315
## 10%         -1.2815516
## 12.5%       -1.1503494
## 15%         -1.0364334
## 17.5%       -0.9345893
## 82.5%        0.9345893
## 85%          1.0364334
## 87.5%        1.1503494
## 90%          1.2815516
## 92.5%        1.4395315
## 95%          1.6448536
## 97.5%        1.9599640
## 99.5%        2.5758293
CuantilesA <- matrix(0,8,2)
CuantilesD <- matrix(0,8,2)
colnames(CuantilesA) <- c('LimInf','LimSup')
colnames(CuantilesD) <- c('LimInf','LimSup')
rownames(CuantilesA) <- c('65','70','75','80','85','90','95','99')
rownames(CuantilesD) <- c('65','70','75','80','85','90','95','99')
CuantilesA[1,1] <-CuantilesData[8]; CuantilesA[1,2] <-CuantilesData[9]
CuantilesA[2,1] <-CuantilesData[7]; CuantilesA[2,2] <-CuantilesData[10]
CuantilesA[3,1] <-CuantilesData[6]; CuantilesA[3,2] <-CuantilesData[11]
CuantilesA[4,1] <-CuantilesData[5]; CuantilesA[4,2] <-CuantilesData[12]
CuantilesA[5,1] <-CuantilesData[4]; CuantilesA[5,2] <-CuantilesData[13]
CuantilesA[6,1] <-CuantilesData[3]; CuantilesA[6,2] <-CuantilesData[14]
CuantilesA[7,1] <-CuantilesData[2]; CuantilesA[7,2] <-CuantilesData[15]
CuantilesA[8,1] <-CuantilesData[1]; CuantilesA[8,2] <-CuantilesData[16]
CuantilesD[1,1] <-Cuantilillos[8];  CuantilesD[1,2] <-Cuantilillos[9]
CuantilesD[2,1] <-Cuantilillos[7];  CuantilesD[2,2] <-Cuantilillos[10]
CuantilesD[3,1] <-Cuantilillos[6];  CuantilesD[3,2] <-Cuantilillos[11]
CuantilesD[4,1] <-Cuantilillos[5];  CuantilesD[4,2] <-Cuantilillos[12]
CuantilesD[5,1] <-Cuantilillos[4];  CuantilesD[5,2] <-Cuantilillos[13]
CuantilesD[6,1] <-Cuantilillos[3];  CuantilesD[6,2] <-Cuantilillos[14]
CuantilesD[7,1] <-Cuantilillos[2];  CuantilesD[7,2] <-Cuantilillos[15]
CuantilesD[8,1] <-Cuantilillos[1];  CuantilesD[8,2] <-Cuantilillos[16]
print(CuantilesD)
##        LimInf    LimSup
## 65 -0.9345893 0.9345893
## 70 -1.0364334 1.0364334
## 75 -1.1503494 1.1503494
## 80 -1.2815516 1.2815516
## 85 -1.4395315 1.4395315
## 90 -1.6448536 1.6448536
## 95 -1.9599640 1.9599640
## 99 -2.5758293 2.5758293
print(CuantilesA)
##        LimInf    LimSup
## 65 -0.8628952 0.8884901
## 70 -0.9319289 1.0207444
## 75 -1.0144701 1.1904247
## 80 -1.1171230 1.3814634
## 85 -1.2529569 1.6102025
## 90 -1.4542681 1.9460444
## 95 -1.8520407 2.4110529
## 99 -1.8520407 3.0342081

Histogramas

col_sequence <- rainbow(n = 7, alpha = 0.35, start = 0, end = 1)
par(mfrow=c(2,1))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2vsumasM3 pEhEx - DATA', lty = 9)

abline(v=CuantilesA[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesA[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesA[2,1], lty=2, col="darkblue"); # 70% INFERIOR
abline(v=CuantilesA[2,2], lty=2, col="darkblue") # 70% SUPERIOR
abline(v=CuantilesA[3,1], lty=2, col="aquamarine4"); # 75% INFERIOR
abline(v=CuantilesA[3,2], lty=2, col="aquamarine4"); # 75% SUPERIOR
abline(v=CuantilesA[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesA[4,2], lty=2, col="green"); # 80% SUPERIOR
abline(v=CuantilesA[5,1], lty=2, col="brown"); # 85% INFERIOR
abline(v=CuantilesA[5,2], lty=2, col="brown"); # 85% SUPERIOR
abline(v=CuantilesA[6,1], lty=2, col="red");  # 90% INFERIOR
abline(v=CuantilesA[6,2], lty=2, col="red"); # 90% SUPERIOR
abline(v=CuantilesA[7,1], lty=2, col="blue");  # 95% INFERIOR
abline(v=CuantilesA[7,2], lty=2, col="blue"); # 95% SUPERIOR
abline(v=CuantilesA[8,1], lty=2, col="orange");  # 99% INFERIOR
abline(v=CuantilesA[8,2], lty=2, col="orange"); # 99% SUPERIOR
legend("topright",
       legend=c("65%","70%","75%","80%","85%","90%","95%","99%"), 
       pch=c(1,2,3,4,5,6,7,8),
       col=c("darkgoldenrod4","darkblue","aquamarine4",
             "green", "brown","red","blue","orange"))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2vsumasM3   pEhEx - ADJUSTED', lty=9)

abline(v=CuantilesD[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesD[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesD[2,1], lty=2, col="darkblue");  # 70% INFERIOR
abline(v=CuantilesD[2,2], lty=2, col="darkblue"); # 70% SUPERIOR
abline(v=CuantilesD[3,1], lty=2, col="aquamarine4");  # 75% INFERIOR
abline(v=CuantilesD[3,2], lty=2, col="aquamarine4"); # 75% SUPERIOR
abline(v=CuantilesD[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesD[4,2], lty=2, col="green"); # 80% SUPERIOR
abline(v=CuantilesD[5,1], lty=2, col="brown");  # 85% INFERIOR
abline(v=CuantilesD[5,2], lty=2, col="brown"); # 85% SUPERIOR
abline(v=CuantilesD[6,1], lty=2, col="red");  # 90% INFERIOR
abline(v=CuantilesD[6,2], lty=2, col="red"); # 90% SUPERIOR
abline(v=CuantilesD[7,1], lty=2, col="blue");  # 95% INFERIOR
abline(v=CuantilesD[7,2], lty=2, col="blue"); # 95% SUPERIOR
abline(v=CuantilesD[8,1], lty=2, col="orange");  # 99% INFERIOR
abline(v=CuantilesD[8,2], lty=2, col="orange"); # 99% SUPERIOR
legend("topright",
       legend=c("65%","70%","75%","80%","85%","90%","95%","99%"), 
       pch=c(1,2,3,4,5,6,7,8),
       col=c("darkgoldenrod4","darkblue","aquamarine4",
             "green", "brown","red","blue","orange"))

Grafica Cuantiles del \(65\%\) y \(80\%\)

par(mfrow=c(2,1))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2vsumasM3  pEhEx - DATA', lty=9)

abline(v=CuantilesA[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesA[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesA[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesA[4,2], lty=2, col="green"); # 80% SUPERIOR
legend("topright",legend=c("65%","80%"), 
       pch=c(1,2),col=c("darkgoldenrod4","green"))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2vsumasM1  - ADJUSTED', lty=9)

abline(v=CuantilesD[1,1], lty=2, col="darkgoldenrod4"); # 65% INFERIOR
abline(v=CuantilesD[1,2], lty=2, col="darkgoldenrod4"); # 65% SUPERIOR
abline(v=CuantilesD[4,1], lty=2, col="green");  # 80% INFERIOR
abline(v=CuantilesD[4,2], lty=2, col="green"); # 80% SUPERIOR
legend("topright",legend=c("65%","80%"), 
       pch=c(1,2),col=c("darkgoldenrod4","green"))

Grafica Cuantiles del \(70\%\) y \(85\%\)

par(mfrow=c(2,1))
hist(tst, breaks = nbreaks, col = col_sequence, 
     main = 'Normalized Log2vsumasM3 pEhEx - DATA', lty=9)

abline(v=CuantilesA[2,1], lty=2, col="darkblue"); 
abline(v=CuantilesA[2,2], lty=2, col="darkblue")
abline(v=CuantilesA[5,1], lty=2, col="brown"); 
abline(v=CuantilesA[5,2], lty=2, col="brown")
legend("topright",legend=c("70%","85%"),
       pch=c(1,2),#3,4,5,6,7,8),
       col=c("brown"))
hist(tst, breaks = nbreaks, col = col_sequence,
     main = 'Normalized Log2vsumasM3  pEhEx - ADJUSTED', lty=9)

abline(v=CuantilesD[2,1], lty=2, col="darkblue"); 
abline(v=CuantilesD[2,2], lty=2, col="darkblue")
abline(v=CuantilesD[5,1], lty=2, col="brown"); 
abline(v=CuantilesD[5,2], lty=2, col="brown")
legend("topright",legend=c("70%","85%"),
       pch=c(1,2),#3,4,5,6,7,8),
       col=c("brown"))